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SUMMARY 

A two-dimensional finite difference computer model for thermal analysis 
of rocket thrust chambers has been developed. The model uses an Iterative 
scheme for calculating the temperature distribution within the chamber wall and 
Implements a successive overrelaxation formula for a quick convergence. The 
Inputs of the model are the dimensions of the thrust chamber wall, types of 
materials used, heat transfer coefficients and temperatures of the hot gas and 
the coolant. The resulting output of the program consists of the nodal temper- 
ature distribution, heat transfer to the coolant and heat transfer from the hot 
gas. 


INTRODUCTION 


The purpose of this report is to describe the two-dimensional finite dif- 
ference model developed for evaluation of the temperature distribution in the 
wall of a rocket thrust chamber. The coating on the inside surface of the 
chamber is subject to high combustion gas temperatures, while the outer walls 
of the nozzle are exposed to space where only radiation heat transfer must be 
considered. A coolant flows through a series of cooling channels located 
inside the chamber wall to maintain the entire component at reasonably moderate 
temperatures. 


The first part of this report describes the finite difference method used 
for the present analysis. Next, the subroutine is described, and, finally 
results of a sample run are presented. 
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NOMENCLATURE 

heat transfer coefficient 

conductivity 

radius 


Ri.R 2 .R 3 .R 4 thermal resistances 


T} , j temperature at node (1,j) 
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radial mesh size 


AR 

A<t> circumferential mesh size 

e relative error for terminating iterations 

eq emisslvity of the outer surface 

a Stefan Boltzman coefficient 

(o successive overrelaxation coefficient 

Subscripts 

A material A 

B material B 

c coolant 

g gas 

1 node 1 

j node j 

o ambient 

Superscripts 

n Iteration number 


THE FINITE DIFFERENCE MODEL 

Consider the geometry of a rocket thrust chamber cross section as shown in 
figure 1. The thrust chamber is made of three materials: nickel, copper, and 

a soot coating. It also consists of a number of cooling channels. Because of 
the symmetry of the configuration, only half of a cooling channel cell is con- 
sidered here (see fig. 2). Since no heat is transferred through the two sides 
of the cell, they are assumed to be Insulated. A finite difference grid is 
superimposed on the aforementioned cell as Is shown in figure 3. The finite 
difference equations Is then written for each nodal point in terms of Its 
neighboring nodes (ref. 1). The finite difference equation for a node located 
in the middle of a material Is given by 
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Similar equations are derived for other nodes, l.e., boundary nodes and nodes 
at the Interface between two different materials. These equations are given In 
appendix A. In general, the finite difference equations give the temperature 
of each node In terms of the temperatures of neighboring nodes and/or heat 
transfer coefficients, conductivities, hot gas, and coolant temperatures. 


The thermal conductivities of the materials used for the combustion wall 
are functions of temperature and these functions are represented by the curves 
In figures 4 to 6 (ref. 3). The gas and coolant temperatures, and the heat 
transfer coefficients for both the hot gas side and the coolant side are also 
known . 


The following numerical procedure Is executed to obtain the temperature 
distribution: 


(1) Set n - 0 

<2> Assume a temperature for each node 

(3) Set n =« n + 1 

(4) Find thermal conductivities based on the nodal temperatures (using 
the data given In figs. 4 to 6) 

(5) Substitute for the nodal temperatures, conductivities and heat trans- 
fer coefficients in the right side of the finite difference equations and 
obtain a new temperature distribution 

(6) If 
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Go to step 8 otherwise, go to step 7 


(7) Use a successive overtaxation formula (ref. 2) to revise the tem- 
perature distribution for a quick convergence. Then go to step 3. 

(8) Stop the calculation, the latest temperature distribution is the 
answer. The successive overrelaxation formula (ref. 2) used in step 7 is given 
by 
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This equation revises the temperature distribution for a quick convergence. 

The most efficient value of to for the geometry under consideration here is 
1.9. This value of <o is obtained by a trial and error procedure to minimize 
the computation time. It should be noted that the successive overrelaxation 
equation makes the convergence four times faster that when it is eliminated 
from the calculation for the configuration considered here. 


THE COMPUTER PROGRAM 

A listing of the computer program developed based on the finite difference 
model discussed in the previous section Is given In the appendix B. The Inputs 
to the program are the dimensions of the combustion chamber, hot gas and cool- 
ant temperatures and heat transfer coefficients, types of materials used at 
each layer, and number of nodes in the radial and circumferential directions. 
Conductivities of five materials are Included in the program. These materials 
are: 

(1) Copper 

(2) Nickel 

(3) Narloy-Z 

(4) Columbium 

(5) Soot (carbon) 


Conductivities of these materials as functions of temperatures are given in 
separate subroutines that are referenced In the main program. 

The output of the program Includes a listing of the temperatures for each 
node, number of Iterations, and the rate of heat transfer from the hot gas to 
the cooling channel. A sample output of the program is Included in appendix B. 

The finite difference program discussed in this report will be used as a 
subroutine in the three dimensional rocket temperature evaluation program. The 
present program provides the rate of heat transfer from hot gas to the coolant 
at each station. Utilizing an iterative scheme and repeating the procedure for 
each station, the temperature distribution and rate of heat transfer to the 
cooling channel can be determined for the whole thrust chamber 
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APPENDIX A 


where 


This appendix presents the finite difference equations (nodal balance of 
energy equations) for the different type nodes. 
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Upper boundary node 


where 



R A4> i 1 


1 


' 3 " AR v k i,j + k i-i,j. 


4 h,R A<t> 
oo 

For forced convection h 0 Is known 


for free convection 


0. 53h 


air 


2R, 


9l T 1,j- T 0 I (2R 0 > Pr 


0.25 


for radiation to space 


'o- V (T U * T o T U * T 1 ,;) T o + T °’ 


6 


Upper boundary (left corner) 
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Right boundary 
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Lower boundary node 


where 
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Lower left side boundary node 
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Lower right side boundary node 



u 


T i J=1 , T i-i.i + Is 

R] R 2 R 3 

i_ 1_ 1_ 

R 1 + R 2 R 3 


where 


AR 


'1 " 2(R + AR/2)A<|> \^k. ^ T k. j _ 1 


R A<t> ( 1 1_ 

AR l k 1,j + k 1-1,j. 


( 3 = h R, A<D 




where 


1 .j 


T I>1„1 T l,3-1 t T I-1J „ 1c 

R, R 2 R 3 R 4 

1 _ 1 _ 1 _ 1 _ 

R 1 + R 2 + R 3 + R 4 


r _ LM. - 
K l ' AR k 


1 1 
— + k 

1.J 


AR 


l 2 * 2 (R + AR/ 2 )A<|> \K ] ^ + k j ^ 


R A 4 > ( 1 ]_\ 

AR yk^ j \_ijJ 

_J 

h c R A 4 » 


14 



Cooling channel side wall 
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Interface between two materials (left boundary node) 



where 


i ,3 


IhlA + lliM + lUzi 

R 1 R 2 R 3 

i_ 1_ !_ 

R 1 + R 2 + R 3 


R_M 


AR A /[l/k A + 1/k. 

*1,3 1 + 1, j. 


+ AR D /| l/k B l/k B 


1,3 1+1 ,3> 


AR, 


1 1 
- + 


'2 2(R AR a / 2)A4» k A T k 


1,3 1 ,3+1- 


AR r 


, B / J _J ' 

( 3 " (R + AR d / 2)A4> k n + k n 

B \ B i J B U-1. 


18 





Interface between two materials (upper left side of the cooling channel) 
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APPENDIX B 


Computer Program Listing and Sample Output 


DATA CASE/'CASE'/,DIME/' SRB9 7 
C 

C FORTRAN PROGRAM FOR EVALUATION OF TEMPERATURE 
C DISTRIBUTION IN SPACECRAFT NOZZLE 
C 

C DEFINITION OF A CELL FOR REGENERATIVE COOLING MODE: 

C A CELL IS THE LEFT HALF OF THE WEDGE CENTERED ON A COOLANT CHANNEL 

C (IE. THE RIGHT SIDE IS THE CENTERLINE OF THE CHANNEL AND THE LEFT 

C SIDE IS THE CENTERLINE OF THE LAND) 

C 

C DEFINITION OF CELL FOR RADIATION COOLING MODE: 

C A CELL IS ONE HALF THE SECTION INTO WHICH THE NOZZLE WALL IS 
C DIVIDED BY THE "NUMBER OF COOLANT CHANNELS" VARIABLE (NCC) 

C 

C I CCW FULL COOLANT CHANNEL WIDTH (FEET) 

C I COOL METHOD OF COOLING 
C =1 REGENERATIVE COOLING 

C =2 RADIATION COOLING 

C CW RADIANS/FULL CELL WIDTH 
C DIFPHIC RADIANS/CELL CHANNEL WIDTH 
C DIFPHIL RADIANS/CELL LAND WIDTH 

C DPHIB CIRCUMFERENTIAL INCREMENT AT CIRCUM INTERFACE(RADIANS) 

C DPHIC CIRCUMFERENTIAL INCREMENT WITHIN CHANNEL (RADIANS) 

C DPHIL CIRCUMFERENTIAL INCREMENT WITHIN LAND (RADIANS) 

C I EM EMMISSIVITY OF CLOSE-OUT MATERIAL 

C I HC1 COOLANT HEAT TRANSFER COEFFICIENT(INNER WALL) 

C I HC2 COOLANT HEAT TRANSFER COEFFICIENT(SIDE WALL) 

C I HC3 COOLANT HEAT TRANSFER COEFFICIENT(OUTER WALL) 

C I HI HOT GAS HEAT TRANSFER COEFFICIENT(BTU/FT**2*S*R) 

C I HOI KNOWN OUTSIDE HEAT TRANSFER COEFFICIENT(BTU/FT**2*S*R) 

C I IFILE FLAG TO MAKE CONTOURING FILE 
C =0 DO NOT CREATE FILE 

C =1 CREATE FILE 

C I IHOUT TYPE OF HEAT TRANSFER AT OUTER BOUNDARY 

C =1 FOR KNOWN OUTSIDE HEAT TRANSFER COEFFICIENT 

C =2 FOR NATURAL CONVECTION AT THE OUTER BOUNDARY 

C =3 FOR RADIATION AT THE OUTER BOUNDARY 

C I MTCH CHANNEL (TOP AND BOTTOM) MATERIAL 

C I MTCLO CLOSE-OUT MATERIAL 
C I MTCOAT COATING MATERIAL 
C I NCC NUMBER OF COOLING CHANNELS 

C I NPHIC NUMBER OF CIRCUM NODES WITHIN CHANNEL AREA 

C I NPHIL NUMBER OF CIRCUM NODES IN LAND AREA (INC. CHANNEL WALL 

C NPHITOT TOTAL NUMBER OF CIRCUMFERENTIAL NODES 

C I NRCHB NUMBER OF RADIAL NODES IN CHANNEL (BOTTOM PORTION) 

C I NRCHT NUMBER OF RADIAL NODES IN CHANNEL (TOP PORTION) 

C I NRCLO NUMBER OF RADIAL NODES IN CLOSE-OUT 

C I NRCOAT NUMBER OF RADIAL NODES IN THE COATING 
C I RCI INNER CHANNEL RADIUS (FEET) 

C I RCO OUTER CHANNEL RADIUS (FEET) 

C I Rl INNER RADIUS (FEET) 

C I RO OUTER RADIUS (FEET) 
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C I TC COOLANT TEMPERATURE (DEG R) 

C I TCOAT COATING THICKNESS (FEET) 

C I Tl HOT GAS TEMPERATURE (DEG R) 

C I TO OUTSIDE TEMPERATURE (DEG R) 

C MATERIAL NUMBERS: (1) COPPER, (2) NICKEL, (3) SOOT (ORIGINAL) 
C (4) NARLOY-Z, (5) RSR 995-AE, (6) COLUMBIUM 
WRITE(6,1)CASE,DIME 
1 FORMAT(1H1///10X,2A4) 

COOL=1 
IHOUT = 2 
IFILE = 0 
NPHIL = 4 
NPHIC = 3 
NRCLO = 5 
NRCHT = 5 
NRCHB =4 
NRCOAT = 0 
MTCLO = 2 

C IF RADIATION COOLING, MTCH MUST BE SAME AS MTCLO 
MTCH = 1 
MTCOAT = 3 
CCW = 0 02/12. 

RO = 1 .945/12. 

RCO = 1.925/12. 

RCI = 1.76/12. 

Rl = 1.74/12. 

TCOAT = 0.000/1 2. 

IF(TCOAT.EQ.O.O)NRCOAT = 0 
NCC = 100 
TO = 530. 

EM =0.9 
Tl = 4748. 

HI =0.01419 
TC = 90. 

HC1 =0.71424 

HC2 = 1.0*0.75406 

HC3 = 1 0*0.77978 

HOI =0.0 

PI = 3.1415 

CW = PI/NCC 

DIFPHIC = CCW/(2.*RCI) 

DIFPHIL = CW-DIFPHIC 
DPHIL = DIFPHIL/(NPHIL-1.) 

IF (NPHIC EQ.99) GOTO 20 
DPHIC = DIFPHIC/NPHIC 
NPHITOT = NPHIC + N PHIL 
GOTO 25 
20 DPHIC = 0. 

NPHITOT = NPHIL 
25 DPHIB = (DPHIC + DPHIL)/2. 

IF (COOL EQ 2) GOTO 15 

WRITE (6,5) N PHIL, NPHIC, CW, DIF PHIL, DIF PHIC.DPHIL, DPHIC, DPHIB 
WRITE(6,2)NRCLO, NRCHT, NRCHB, NRCOAT 
GOTO 10 
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15 WRITE (6,8)NPHIL,CW,DPHIL 
WRITE (6,11)MTCL0,NRCL0 
WRITE (6,13)MTCOAT,NRCOAT 

10 CONTINUE 
C 

2 FORMAT(///10X,'NRCLO = ',I5,2X,'NRCHT = ',I5,2X,'NRCHB = ',I5,2X, 

1'NRC0AT = ',I5) 

5 FORMAT (//10X,'NPHIL = ',I4,3X,'NPHIC = ',I4,3X./10X,'CW = ',F8.6. 

13X,'DIFPHIL = ',F8.6,3X,'DIFPHIC = ',F8.6,3X,/10X,'DPHIL = ',F8.6, 

23X/DPHIC = ',F8.6,3X,'DPHIB = ',F8.6) 

8 FORMAT (//10X/NPHI = ',14,/IOX/CELL WIDTH = ',F8 6,/10X,'DPHI = 

1F8.6) 

11 FORMAT (//10X, 'SUBSTRATE: MATERIAL = ',I3,6X,'N0DES = ',I3) 

13 FORMAT (/10X, 'COATING: MATERIAL^ ',13, 6X, NODES = ',I3) 

C 

CALL COND(RI,RO,DPHIC,DPHIL,DPHIB,NRCLO,NRCHT,NRCHB,NRCOAT, 
1NPHIC,NPHIL,NPHITOT,RCI,RCO,TCOAT,TO,EM,TI,HI,TC,HC1,HC2,HC3, 
2MTCL0,MTCH,MTC0AT,C00L,IH0UT,H01,IFILE,NCC) 

STOP 

END 

C 

SUBROUTINE COND(RI,RO,DPHIC,DPHIL,DPHIB,NRCLO,NRCHT,NRCHB,NRCOAT, 
1NPHIC,NPHIL,NPHITOT,RCI,RCO,TCOAT,TO,EM,TI,HI,TC,HC1,HC2,HC3, 

2MTCLO, MTCH.MTCOAT, COOL, IHOUT, HOI, IFILE.NCC) 

C 

C COOL METHOD OF COOLING (FROM $MAIN) 

C Cl -2-3 (L OR C) GIVES GEOMETRY TO RESISTANCES IN AREA/LENGTH 
C DOT ARRAY WITH SCALING DATA FOR CONTOURING ROUTINE 
C DPHI CIRCUMFERENTIAL INCREMENT (GIVEN C OR L VALUE) 

C DRCHB RADIAL INCREMENT IN BOTTOM OF CHANNEL (FROM $MAIN) 

C DRCHT RADIAL INCREMENT IN TOP OF CHANNEL (FROM $MAIN) 

C DRCLO RADIAL INCREMENT IN CLOSE-OUT (FROM $MAIN) 

C DRCOAT RADIAL INCREMENT IN COATING (FROM $MAIN) 

C EE FRACTIONAL ERROR AT END OF CONDUCTIVITY ITERATION 

C EM EMMISSIVITY OF CLOSE-OUT MATERIAL (FROM SMAIN) 

C ERR FRACTIONAL ERROR WITHIN CONDUCTIVITY ITERATION 

C ERROR RELATIVE ERROR FOR TERMINATING ITERATIONS 

C FC() ARRAY WITH CONTOURING LEVELS FOR CONTOURING ROUTINE 
C GBROV = G'BETA‘ROU"2/VISC**2 

C HC1-2-3 COOLANT HEAT TRANSFER COEFFICIENT (FROM SMAIN) 

C HI INSIDE HEAT TRANSFER COEFFICIENT (FROM SMAIN) 

C HO( ) OUTSIDE HEAT TRANSFER COEFFICIENT FOR EACH NODE 
C HOAVE AVERAGE OUTSIDE HEAT TRANSFER COEFFICIENT 
C HOI KNOWN OUTSIDE HEAT TRANSFER COEFFICIENT (FROM SMAIN) 

C ICH( , ) HOLDS TYPE OF POINT 

C IDASH ARRAY WITH DATA FOR CONTOURING ROUTINE 
C IEE ERROR COUNTER AT END OF CONDUCTIVITY ITERATION 
C IERR ERROR COUNTER WITHING CONDUCTIVITY ITERATION 
C (FILE FLAG TO MAKE CONTOURING FILE 
C IFLAG() ARRAY WITH DATA FOR CONTOURING ROUTINE 
C IHOUT TYPE OF HEAT TRANSFER AT OUTER BOUNDARY (FROM SMAIN) 

C ITER COUNTER WITHIN SINGLE CONDUCTIVITY ITERATION 
C ITER1 COUNTER: NUMBER OF CONDUCTIVITY ITERATIONS 
C ITER2 COUNTER FOR TOTAL NUMBER OF ITERATIONS 


23 


C LOOP VARS: I,I1,I2,I3,II,J,J1,J2,J3,JJ,K,M 

C MTCH CHANNEL (TOP AND BOTTOM) MATERIAL (FROM SMAIN) 

C MTCLO CLOSE-OUT MATERIAL (FROM SMAIN) 

C MTCOAT COATING MATERIAL (FROM SMAIN) 

C NCC NUMBER OF COOLING CHANNELS 

C NPHIC-L-TOT NUMBER OF CIRCUMFERENTIAL NODES (FROM SMAIN) 

C NCON NUMBER OF CONTOURING LEVELS IN FC() 

C NR NUMBER OF RADIAL NODES (NRCLO + NRCHT + NRCHB + NRCOAT) 
C NRCHB-CHT-CLO-COAT NUMBER RADIAL NODES (FROM SMAIN) 

C NR1 NUMBER OF RADIAL NODES (NRCLO) 

C NR2 NUMBER OF RADIAL NODES (NRCLO + NRCHT) 

C NR3 NUMBER OF RADIAL NODES (NRCLO + NRCHT + NRCHB) 

C OMEGA SOR COEFFICIENT 

C PR PRANDTL NUMBER 

C Q SUM OF ALL THE HEAT FLOWS 

C QC TOTAL HEAT FLOW TO SINGLE COOLANT CHANNEL 

C QC1 HEAT FLOW TO INNER CHANNEL WALL 

C QC2 HEAT FLOW TO SIDE CHANNEL WALL 

C QC3 HEAT FLOW TO UPPER CHANNEL WALL 

C Ql HEAT FLOW THROUGH HOT GAS WALL 

C QO HEAT FLOW TO OUTSIDE 

C RADP( , ) HOLDS RADIUS OF NODES (INCHES) 

C RCI-CO-I-O WALL RADII (FROM SMAIN) 

C RKAIR CONDUCTIVITY OF AIR 

C RK1-2-3-4-5 DUMMY CONDUCTIVITIES (ALSO WITH CH, CLO AND COAT) 

C R1-2-3-4-5 DUMMY THERMAL RESISTANCES 
C SIGMA STEFAN-BOLTZMANN CONSTANT 

C T( , ,1) TEMPERATURE DISTRIBUTION USED FOR CONDUCTIVITIES 
C T( , ,2) CONVERGING TEMPERATURE DISTRIBUTION AT SINGLE COND 
C TC COOLANT TEMPERATURE (FROM SMAIN) 

C TCOAT COATING THICKNESS (FROM SMAIN) 

C Tl HOT GAS TEMPERATURE (FROM SMAIN) 

C TO OUTSIDE TEMPERATURE (FROM SMAIN) 

C TR1-2-3-4-5 "DUMMY "TEMPERATURES TO DETERMINE CONDUCTIVITIES 
C TT2 CONVERGENCE CHECK AT END OF CONDUCTIVITY 
C 
C 

DIMENSION T(40,40,2),ICH(40,40),RADP(40,40),HO(40),IFLAG(7) 

DIMENSION DOT(8),FC(20),IDASH(10) 

C IF NCON LE. 0 THEN ALL PRINTED OUTPUT IS SUPPRESSED 
NCON =-15 
FC(1) = 1300. 

FC(2) = 1250. 

FC(3) = 1200. 

FC(4) = 1150. 

FC(5) = 1100. 

FC(6) = 1050. 

FC(7) = 1000. 

FC(8) = 950. 

FC(9) =900. 

FC(10) = 850. 

FC(1 1) = 800. 

FC(12) = 750. 

FC(13) = 700. 
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FC(14) = 660. 

FC(15) =650. 

IFLAG(I) = NRCLO + NRCHT + NRCHB + NRCOAT 

IFLAG(2) = NPHITOT 

IFLAG(3) = NCON 

IFLAG(4) = 1 

IFLAG(5) = 1 

IFLAG(6) = 0 

IFLAG(7) = 0 

DOT(1) =0. 

DOT(2) =0. 

DOT(3)=0. 

DOT(4) =0. 

DOT(5) =0. 

DOT(6) =0. 

DOT(7) =0. 

DOT(8) =0. 

IDASH(I) = 1 
OMEGA= 1 9 

SIGMA =0.173*10. **{-8)/3600. 

RKAIR = 0.014/3600. 

GBROV = 3. 160*10**6 
PR =0.72 
ERROR = 10.**(-7) 

DRCLO = (RO-RCO)/{NRCLO-1 .) 

IF (COOL .EQ. 1) GOTO 516 
DRCHT = 0. 

DRCHB = 0. 

GOTO 517 

516 DRCHT = (RCO-RCI)/NRCHT 
DRCHB = (RCI-RI)/NRCHB 

517 IF(NRCOAT.EQ.O)GO TO 512 
DRCOAT = TCOAT/NRCOAT 

512 NR1 = NRCLO 

NR2 = NRCLO + NRCHT 

NR3 = NRCLO + NRCHT + NRCHB 

NR = NRCLO + NRCHT + NRCHB + NRCOAT 

IF (COOL EQ 2) GOTO 6 

WRITE(6,2)DRCLO, DRCHT, DRCHB, DRCOAT.DPHIC, DPHIL, NR, NR1,NR2,NR3 
GOTO 7 

6 WRITE (6,9)DRCLO,DRCOAT,NPHIL,NR1 ,NR 

7 CONTINUE 

2 FORMAT(/10X, DRCLO = ',F8. 6, /I OX, 'DRCHT = ',F8.6,/10X,'DRCHB = ', 
1F8.6./10X, 'DRCOAT = ',F8.6,/10X,'DPHIC = ',F8.6,/10X, 'DPHIL = ', 
2F8.6,/10X,'NR = M5,2X,/10X,'NR1=',I5,2X,'NR2 = ',I5,2X,'NR3 = ',I5) 

9 FORMAT (//10X/DRSUBSTRATE = ',F8.6,/10X,'DRCOAT = ',F8.6,/10X, 

1'NPHI = ',I4,//10X,'NR1 =',I4,5X,'NR = ',I4) 

DO 500 II =1, NPHITOT 
HO(II) = HOI 
DO 501 J1 =1,NR 
T(I1,J1,2) =(TI + TO + TC)/3. 

501 CONTINUE 
500 CONTINUE 
ITER1 =0 


ITER2 =0 

4 DO 503 12 = I.NPHITOT 
DO 504 J2 = 1,NR 
T(I2,J2,1) = T(I2,J2,2) 

504 CONTINUE 
503 CONTINUE 

ITER =0 

ITER1 =ITER1 +1 
1 IERR = 0 
ITER = ITER + 1 
ITER2 = ITER2 + 1 
DO 1000 I = I.NPHITOT 
DO 1001 J = 1,NR 
C WRITE(6, 507)1, J 

507 FORMAT(10X,2(I3,2X)) 

C 

IF(J.LT.NR1)GO TO 508 
IF(J.GT.NR1.AND.J.LT.NR2)GO TO 509 
IF(J.GT.NR2.AND.J.LT.NR3)GO TO 510 
IF(J.GT.NR3.AND.TCOAT.NE.O.O)GO TO 505 
IF(J.EQ.NR3.AND.TCOAT.EQ.O.O)GO TO 510 
GO TO 506 

505 R = RI-(J-NR3)*DRCOAT 
DR = DRCOAT 

GO TO 511 

508 R = RO-(J-1 )’DRCLO 
DR = DRCLO 

GO TO 511 

509 R = RCO-(J-NR1)*DRCHT 
DR = DRCHT 

GO TO 511 

510 IF (COOL .EQ. 2) GOTO 508 
R = RCI-(J-NR2)’DRCHB 

DR = DRCHB 

511 C1L = DR/(R'DPHIL) 

C2L = (R-DR/2.)*DPHIL/DR 
C3L = (R + DR/2.)*DPHIL/DR 
IF (COOL EQ. 2) GOTO 515 
C1C = DR/(R*DPHIC) 

C2C = (R-DR/2.)*DPHIC/DR 
C3C = (R + DR/2.)*DPHIC/DR 
C 

515 IF (I.GT.NPHIL) GOTO 513 
Cl =C1L 
C2 ~C2L 
C3 = C3L 
DPHI = DPHIL 
GOTO 506 
513 Cl =C1C 
C2 — C2C 
C3 = C3C 
DPHI = DPHIC 

506 TT2 = T(l, J,2) 

C 
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IF (COOL .EQ. 2) GOTO 400 
IF(I.EQ.NPHIL)GOTO 333 
IF(I.EQ.1.AND.J.EQ.1)GO TO 100 
IF(I.EQ.NPHITOT.AND.J.EQ.1)GO TO 110 
IF(J.EQ.1)GO TO 120 
IF(I.EQ.1.AND.J.EQ.NR)GO TO 130 
IF(I.EQ.NPHITOT.AND.J.EQ.NR)GO TO 140 
IF(J.EQ.NR)GO TO 150 
IF(I.EQ.1.AND.J.EQ.NR3)GO TO 135 
IF(I EQ.NPHITOT.AND.J.EQ.NR3)GO TO 145 
IF(J.EQ.NR3)GO TO 155 
IF(I.EQ.1)GO TO 160 

IF(I.GT.NPHIL.AND.J.EQ.NR1)GO TO 180 
IF(J.GT.NR1.AND.J.LT.NR2)GO TO 190 
GO TO 191 
C 

190 IF(I.GT.NPHIL)GO TO 210 

191 IF(I.GT.NPHIL.AND.J.EQ.NR2)GO TO 230 
IF(J.EQ.NR1)GO TO 240 
IF(J.EQ.NR2)GO TO 260 
IF(I.EQ.NPHITOT)GO TO 250 

GOTO 192 

333 IF (J EQ. 1) GOTO 300 
IF (J EQ. NR1) GOTO 170 
IF (J .EQ. NR2) GOTO 220 
IF (J EQ. NR) GOTO 320 
IF (J EQ. NR3) GOTO 310 
IF ((J .GT. NR1) .AND. (J .LT. NR2)) GOTO 200 

MIDDLE NODE: CIRCUMFERENTIAL INTERFACE 


TR1 =T(I,J,1) 

TR2 = T(I + 1,J,1) 

TR3 = T(I,J-1,1) 

TR4 = T(I-1,J,1) 

TR5 = T(I,J + 1,1) 

IF (J GT. NR3) GOTO 334 
IF (J .GT. NR1) GOTO 335 
CALL METAL (MTCLO, RK1, TR1) 
CALL METAL (MTCLO, RK2, TR2) 
CALL METAL (MTCLO, RK3, TR3) 
CALL METAL (MTCLO, RK4, TR4) 
CALL METAL (MTCLO, RK5, TR5) 
GOTO 336 

334 CALL METAL(MTCOAT, RK1, TR1) 
CALL METAL(MTCOAT, RK2, TR2) 
CALL METAL(MTCOAT, RK3, TR3) 
CALL METAL(MTCOAT, RK4, TR4) 
CALL METAL(MTCOAT, RK5, TR5) 
GOTO 336 

335 CALL METAL (MTCH, RK1, TR1) 
CALL METAL (MTCH, RK2, TR2) 
CALL METAL (MTCH, RK3, TR3) 
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CALL METAL (MTCH, RK4, TR4) 

CALL METAL (MTCH, RK5, TR5) 

336 CONTINUE 

R1 = (1./RK1 + 1./RK2)/(2.*C1C) 

R2 =(1./RK1 + 1./RK3)’DR/(2."DPHIB*(R + DR/2.)) 

R3 = (1./RK1 + 1./RK4)/(2.*C1L) 

R4 = (1./RK1 + 1 /RK5)'DR/(2.‘DPHIB*(R-DR/2.)) 

T(I,J,2) = ((T(I + 1,J,2)/R1 + T(I,J-1,2)/R2 + T(I-1,J,2)/R3 + 

1 T(l, J + 1 ,2)/R4)/(1./R1 + 1 ,/R2 + 1 ,/R3 + 1 ,/R4)-T(l,J,2))*OMEGA + 
2T(I,J,2) 

ICH(I.J) —50 
GO TO 1002 

RADIATION COOLED DECISION SECTION 

400 IF ((I.EQ.I).AND.(J.EQ.I)) GOTO 100 

IF ((I.EQ.NPHITOT).AND.(J.EQ.I)) GOTO 110 
IF (J.EQ.1) GOTO 120 
IF ((I.EQ.1). AND. (J.EQ.NR)) GOTO 130 
IF ((I.EQ.NPHITOT). AND (J.EQ.NR)) GOTO 140 
IF (J.EQ.NR) GOTO 150 
IF ((I.EQ.1).AND.(J.EQ.NR3)) GOTO 135 
IF ((I.EQ.NPHITOT).AND.(J.EQ.NR3)) GOTO 145 
IF (J.EQ.NR3) GOTO 155 
IF (I.EQ.1) GOTO 160 
IF (I.EQ.NPHITOT) GOTO 250 

MIDDLE NODE 

192 TR1 =T(I,J,1) 

TR2=T(I + 1,J,1) 

TR3 =T(I,J-1,1) 

TR4 = T(M,J,1) 

TR5 = T(I.J + 1,1) 

IF(J.GT.NR1.AND.J.LT.NR3)GO TO 101 

IF(J.GT.NR3)GO TO 109 

CALL METAL (MTCLO, RK1.TR1) 

CALL METAL (MTCLO, RK2.TR2) 

CALL METAL (MTCLO, RK3.TR3) 

CALL METAL (MTCLO, RK4.TR4) 

CALL METAL (MTCLO, RK5.TR5) 

GO TO 102 

101 CALL METAL (MTCH, RK1.TR1) 

CALL METAL (MTCH, RK2.TR2) 

CALL METAL (MTCH, RK3.TR3) 

CALL METAL (MTCH, RK4.TR4) 

CALL METAL (MTCH, RK5.TR5) 

GO TO 102 

109 CALL METAL (MTCOAT, RK1.TR1) 

CALL METAL (MTCOAT, RK2.TR2) 

CALL METAL (MTCOAT. RK3.TR3) 

CALL METAL (MTCOAT, RK4.TR4) 

CALL METAL (MTCOAT, RK5.TR5) 

102 R1 = 1 ,/2./C1*(1./RK1 + 1./RK2) 
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R2 = 1./2./C3*(1./RK1 +1./RK3) 

R3 = 1./2./C1*(1./RK1 +1./RK4) 

R4 = 1./2./C2*(1./RK1 +1./RK5) 

T(I,J,2) = ((T(I + 1,J,2)/R1 + T(I,J-1,2)/R2 + T(I-1,J,2)/R3 + 

1T(I,J + 1 ,2)/R4)/(1 ,/R1 + 1 ,/R2 + 1 ,/R3 + 1 ./R4)-T(l,J,2))*OMEGA + 
2T(I,J,2) 

ICH(I,J) = 0 
GO TO 1002 

UPPER LEFT SIDE CORNER BOUNDARY NODE 
100 IF(IHOUT.EQ.1)HO(l) = HOI 

IF(IHOUT EQ.2)HO(l) =0 53*RKAIR/(2.*RO)’(GBROV*ABS(T(l,J,1)-TO) 
1 “(2.'RO)**3'PR)"0.25 

IF(IHOUT.EQ.3)HO(l)=EM*SIGMA'(T(l,J,1)”3.+TO‘T(l,J,1)**2.+ 

1T(l,J,1)*TO**2.+TO**3.) 

TR1 =T(I,J,1) 

TR2 = T(I + 1,J,1) 

TR3 = T(I,J + 1,1) 

CALL METAL (MTCLO, RK1.TR1) 

CALL METAL (MTCLO, RK2.TR2) 

CALL METAL (MTCLO, RK3.TR3) 

R1 = 1./C1*(1./RK1 + 1./RK2) 

R2 = 1./C2*(1./RK1 +1./RK3) 

IF(HO(I).EQ.O.)GO TO 103 
R3 = 2./HO(l)/RO/DPHI 
GO TO 104 

103 R3 = 10. **16 

104 T(I,J,2) = ((T(I + 1,J.2)/R1 +T(I,J + 1.2J/R2 +TO/R3)/ 

1(1./R1 + 1./R2 + 1./R3)-T(l,J,2))'OMEGA +T(I,J,2) 

ICH(I,J) = 100 

GO TO 1002 

UPPER RIGHT HAND SIDE CORNER BOUNDARY NODE 
110 IF(IHOUT.EQ.1)HO(l) = HOI 

IF(IHOUT.EQ.2)HO(I)=0.53*RKAIR/(2.‘RO)*(GBROV'ABS(T(I,J,1)-TO) 

1*(2,*RO)**3*PR)**0.25 

IF(IHOUT,EQ.3)HO(l) = EM*SIGMA*(T(I,J,1)**3. +TO*T(l,J,1)**2. + 
1T(l,J,1)*TO**2. + TO**3.) 

TR1 =T(I,J,1) 

TR2 = T(I-1,J,1) 

TR3 = T(I,J + 1,1) 

CALL METAL (MTCLO, RK1.TR1) 

CALL METAL (MTCLO, RK2.TR2) 

CALL METAL (MTCLO, RK3.TR3) 

R1 = 1./C1*(1./RK1 + 1./RK2) 

R2 = 1./C2*(1./RK1 +1./RK3) 

IF(HO(I).EQ.O.)GO TO 105 
R3 = 2./HO(l)/RO/DPHI 
GO TO 106 

105 R3 = 10.”16 

106 T(I,J,2) = ((T(I-1,J,2)/R1 +T(I,J + 1,2)/R2 + TO/R3)/ 

1(1. /R1 + 1./R2 + 1./R3)-T(l,J,2))-OMEGA +T(I,J,2) 
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ICH(I,J) = 110 
GO TO 1002 

UPPER BOUNDARY NODE 
120 IF(IHOUT EQ.I)HO(I) = H01 

IF(IHOUT.EQ.2)HO(I)=0.53*RKAIR/(2.*RO)*(GBROV*ABS(T(I,J,1)-TO) 
1*(2 'RO)**3'PR)*‘0.25 

IF(IHOUT.EQ.3)HO(l) = EM*SIGMA*(T(l,J,1)**3.+TO’T(l,J,1)**2.+ 
1T(l,J,1)*TO**2.+TO**3.) 

TR1 = T(I,J,1) 

TR2 = T(I + 1,J,1) 

TR3 = T(I,J + 1,1) 

TR4 = T(I-1,J,1) 

CALL METAL (MTCLO, RK1.TR1) 

CALL METAL (MTCLO, RK2.TR2) 

CALL METAL (MTCLO, RK3.TR3) 

CALL METAL (MTCLO, RK4.TR4) 

R1 = 1./C1"(1./RK1 + 1./RK2) 

R2 = 1 ,/2./C2*(1 /RK1 + 1 ./RK3) 

R3 = 1./C1*(1./RK1 + 1./RK4) 

IF(HO(I).EQ.O.O)GO TO 107 
R4 = 1 /HO(!)/RO/DPHI 
GO TO 108 

107 R4 = 10. **16 

108 T(I,J,2) = ((T(l + 1,J,2)/R1 +T(I,J + 1,2)/R2 +T(I-1,J,2)/R3 +TO/R4) 

1/( 1 /R1+ 1 ,/R2 + 1 ,/R3 + 1 ,/R4)-T(l,J,2))*OMEGA + T(I,J,2) 

ICH(I.J) = 120 
GO TO 1002 

UPPER BOUNDARY: CIRCUMFERENTIAL INTERFACE 

300 IF(IHOUT.EQ.1)HO(l) = HOI 

IF(IHOUT EQ.2)HO(I)=0 53*RKAIR/(2.*RO)*(GBROV*ABS(T(l,J,1)-TO) 
1 *(2.'RO)‘“3‘PR)’*0.25 

IF(IHOUT.EQ 3)HO(l) ~EM'SIGMA*(T(I,J,1)'*3. +TO*T(l,J,1)**2 + 
1T(l,J,1)‘TO'*2.+TO**3.) 

TR1 =T(I,J,1) 

TR2 = T(I + 1,J,1) 

TR3 =T(I,J + 1,1) 

TR4 = T(I-1,J,1) 

CALL METAL (MTCLO, RK1.TR1) 

CALL METAL (MTCLO, RK2.TR2) 

CALL METAL (MTCLO, RK3.TR3) 

CALL METAL (MTCLO, RK4.TR4) 

R1=(1 /RK1 + 1./RK2)/C1C 

R2 = (1 /RK1 +1 /RK3)*DRCLO/(2.*DPHIB*(RO-DRCLO/2.)) 

R3 = (1 /RK1 + 1./RK4)/C1L 
IF(HO(l) EQ 0.0)GO TO 301 
R4 = 1./(HO(l)*RO*DPHIB) 

GO TO 302 

301 R4 = 10**16 

302 T(I,J,2) = ((T(l + 1,J,2)/R1 + T(I,J + 1,2)/R2 +T(I-1,J,2)/R3 +TO/R4) 

1/(1 ,/R1 + 1./R2 + 1./R3 + 1./R4)-T(l,J,2))*OMEGA + T(l,J,2) 
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ICH(I,J) = 300 
GO TO 1002 

LOWER LEFT SIDE BOUNDARY NODE 

130 TR1 =T(I,J,1) 

TR2 = T(I + 1,J,1) 

TR3 = T(I,J-1,1) 

IF(TCOAT.EQ.O.O)GO TO 131 
CALL METAL (MTCOAT, RK1.TR1) 

CALL METAL (MTCOAT, RK2.TR2) 

CALL METAL (MTCOAT, RK3.TR3) 

GO TO 132 

131 CALL METAL (MTCH, RK1.TR1) 

CALL METAL (MTCH, RK2.TR2) 

CALL METAL (MTCH, RK3.TR3) 

132 R1 = 1./C1*(1./RK1 + 1./RK2) 

R2 = 1./C3*(1./RK1 + 1./RK3) 

R3 = 2./HI/R/DPHI 

T( I , J ,2) =((T(I + 1,J,2)/R1 + T(I,J-1,2)/R2 + TI/R3)/ 

1(1 ,/R 1 + 1./R2 + 1,/R3)-T(l,J,2))*OMEGA + T(I,J,2) 

ICH(I,J) = 130 
GO TO 1002 

INTERFACE BETWEEN COPPER AND COATING(LEFT BOUNDARY) 

135 R = RI 
TR1 =T(I,J,1) 

TR2 = T(I + 1,J,1) 

TR3 = T(I,J + 1,1) 

TR5 = T(I,J-1,1) 

CALL METAL (MTCH, RKCH1.TR1) 

CALL METAL (MTCOAT, RKCOA1.TR1) 

CALL METAL (MTCH, RKCH2.TR2) 

CALL METAL (MTCOAT, RKCOA2.TR2) 

CALL METAL (MTCOAT, RKCOA3.TR3) 

CALL METAL (MTCH, RKCH5.TR5) 

R1 = R*DPHI/(DRCHB/(1 ./RKCH1 + 1 ./RKCH2) + DRCOAT/(1./RKCOA1 + 1./RKCOA2)) 
R2 = DRCOAT/1./(R-DRCOAT/2.)/DPHI*(1./RKCOA1 + 1./RKCOA3) 

R4 = DRCHB/1./(R + DRCHB/2.)/DPHI*(1./RKCH1 +1./RKCH5) 

T(I,J,2) =((1.*T(I + 1,J,2)/R1+T(I,J + 1,2)/R2 + T(I,J-1,2)/R4 
1 )/( 1 ,/R1 + 1 /R2 + 1./R4)-T(l,J,2))*OMEGA + T(I,J,2) 

ICH(I,J) = 135 
GO TO 1002 

LOWER RIGHT SIDE BOUNDARY NODE 

140 TR1 =T(I,J,1) 

TR2 = T(I-1,J,1) 

TR3 = T(I,J-1,1) 

IF(TCOAT.EQ.O.O)GO TO 141 
CALL METAL (MTCOAT, RK1.TR1) 

CALL METAL (MTCOAT, RK2.TR2) 

CALL METAL (MTCOAT, RK3.TR3) 
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GO TO 142 

141 CALL METAL (MTCH, RK1.TR1) 

CALL METAL (MTCH, RK2.TR2) 

CALL METAL (MTCH, RK3.TR3) 

142 R1 = 1./C1*(1,/RK1 + 1./RK2) 

R2 = 1./C3*(1./RK1 + 1./RK3) 

R3 = 2./HI/R/DPHI 

T(I,J,2) = ((T(I-1,J,2)/R1 + T(I,J-1 ,2)/R2 + TI/R3)/ 

1(1./R1 + 1./R2 + 1./R3)-T(l,J,2))*OMEGA + T(I,J,2) 

ICH(I,J) = 140 
GO TO 1002 

INTERFACE BETWEEN COPPER AND COATING(RIGHT BOUNDARY) 

145 R = RI 
TR1 = T(I,J,1) 

TR3 = T(I,J + 1.1) 

TR4 = T(I-1 ,J,1) 

TR5 = T(I,J-1,1) 

CALL METAL (MTCH, RKCH1.TR1) 

CALL METAL (MTCOAT, RKCOA1.TR1) 

CALL METAL (MTCOAT, RKCOA3.TR3) 

CALL METAL (MTCH, RKCH4.TR4) 

CALL METAL (MTCOAT, RKCOA4.TR4) 

CALL METAL (MTCH, RKCH5.TR5) 

R2 = DRCO AT/1. /(R-DRCO AT/2. )/DPHI*(1./RKCOA1 + 1./RKCOA3) 

R3 = R'DPHI/(DRCHB/(1 ./RKCH1 + 1 ./RKCH4) + DRCOAT/(1 ,/RKCOAI + 1 ./RKCOA4)) 
R4 = DRCHB/1./(R + DRCHB/2.)/DPHI*(1./RKCH1 +1./RKCH5) 

T(I,J,2) =((T(I,J + 1,2)/R2 + 1 *T(I-1,J,2)/R3 + T(I,J-1,2)/R4 
1)/(1./R2 + 1./R3 + 1./R4)-T(l,J,2))*OMEGA + T(l,J,2) 

ICH(I.J) = 145 
GO TO 1002 

LOWER BOUNDARY NODES 

150 TR1 =T(I,J,1) 

TR2 = T(I + 1 , J, 1) 

TR3 = T(I,J-1,1) 

TR4 = T(I-1,J,1) 

IF(TCOAT.EQ.O.O)GO TO 151 
CALL METAL (MTCOAT, RK1.TR1) 

CALL METAL (MTCOAT, RK2.TR2) 

CALL METAL (MTCOAT, RK3.TR3) 

CALL METAL (MTCOAT, RK4.TR4) 

GO TO 152 

151 CALL METAL (MTCH, RK1.TR1) 

CALL METAL (MTCH, RK2.TR2) 

CALL METAL (MTCH, RK3.TR3) 

CALL METAL (MTCH, RK4.TR4) 

152 R1 = 1 /C1*(1./RK1 + 1./RK2) 

R2 — 1./2./C3*(1,/RK1 + 1./RK3) 

R3 = 1./C1*(1./RK1 + 1./RK4) 

R4 = 1./HI/R/DPHI 

T(I,J,2) = ((T(I + 1 , J ,2)/R 1 +T(I,J-1,2)/R2 + T(I-1,J,2)/R3 +TI/R4)/ 
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1< 1 ,/R1 + 1 ./ R2 + 1 ,/R3 + 1 ,/R4)-T(l,J,2))*OMEGA + T(l, J,2) 

ICH(I,J) = 150 
GO TO 1002 

LOWER BOUNDARY NODE: CIRCUMFERENTIAL INTERFACE 

10 TR1 =T(I,J,1) 

TR2 = T(I + 1,J,1) 

TR3 = T(I,J-1,1) 

TR4 = T{1-1,J,1) 

IF(TCOAT.EQ.O.O)GO TO 321 
CALL METAL (MTCOAT, RK1.TR1) 

CALL METAL (MTCOAT, RK2.TR2) 

CALL METAL (MTCOAT, RK3.TR3) 

CALL METAL (MTCOAT, RK4.TR4) 

GO TO 322 

!l CALL METAL (MTCH, RK1.TR1) 

CALL METAL (MTCH, RK2.TR2) 

CALL METAL (MTCH, RK3.TR3) 

CALL METAL (MTCH, RK4.TR4) 

2 R1 =(1./RK1 + 1./RK2)/C1C 
R2 = (1./RK1 4- 1 ./RK3)*DR/(2.*DPHIB*(R + DR/2.)) 

R3 = (1./RK1 + 1./RK4)/C1L 
R4 = 1 /(HI'R'DPHIB) 

T(I,J,2) = ((T(l + 1 ,J,2)/R1 + T(l, J-1 ,2)/R2 + T(l-1 , J,2)/R3 + TI/R4)/ 
1(1 ,/R1 + 1 ,/R2 + 1 ,/R3 + 1./R4)-T(l,J,2))’OMEGA + T(I,J,2) 

ICH(«, J) — 320 
GO TO 1002 


INTERFACE BETWEEN COPPER AND COATING 

155 R = RI 
TR1 =T(I,J,1) 

TR2 — T(I + 1,J,1) 

TR3 =T(I,J + 1,1) 

TR4 = T(I-1,J,1) 

TR5 = T(I,J-1,1) 

CALL METAL (MTCH, RKCH1.TR1) 

CALL METAL (MTCOAT, RKCOA1.TR1) 

CALL METAL (MTCH, RKCH2.TR2) 

CALL METAL (MTCOAT, RKCOA2.TR2) 

CALL METAL (MTCOAT, RKCOA3.TR3) 

CALL METAL (MTCH, RKCH4.TR4) 

CALL METAL (MTCOAT, RKCOA4.TR4) 

CALL METAL (MTCH, RKCH5.TR5) 

R1 = R*DPHI/(DRCHB/(1 ./RKCH1 + 1 ./RKCH2) + DRCOAT/(1 ,/RKCOAI + 1 ./RKCOA2)) 
R2 = DRCOAT/2 /(R-DRCOAT/2.)/DPHI*(1 ,/RKCOAI + 1 ./RKCOA3) 

R3 = R“DPHI/(DRCHB/(1 ,/RKCHI + 1 ./RKCH4) + DRCOAT/(1 ,/RKCOAI + 1 ./RKCOA4)) 
R4 - DRCHB/2./(R + DRCHB/2.)/DPHI*(1 ,/RKCHI + 1 /RKCH5) 

T(I,J,2) - ((T(l + 1,J,2)/R1 +T(I,J + 1,2)/R2 + T(I-1,J,2)/R3 + T(I,J-1,2) 

1/R4)/(1./R1 + 1 ,/R2 + 1 ,/R3 + 1./R4)-T(l,J,2))*OMEGA + T(I,J,2) 

ICH(I.J) = 155 
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GO TO 1002 


CHANNEL-COATING INTERFACE: CIRCUMFERENTIAL INTERFACE 

310 R = Rl 
TR1 — T(I,J,1) 

TR2 = T(I-M,J,1) 

TR3 =T(I,J + 1,1) 

TR4 =T(I-1,J,1) 

TR5 = T(I,J-1,1) 

CALL METAL (MTCH, RKCH1.TR1) 

CALL METAL (MTCOAT, RKCOA1.TR1) 

CALL METAL (MTCH, RKCH2.TR2) 

CALL METAL (MTCOAT, RKCOA2.TR2) 

CALL METAL (MTCOAT, RKCOA3.TR3) 

CALL METAL (MTCH, RKCH4.TR4) 

CALL METAL (MTCOAT, RKCOA4.TR4) 

CALL METAL (MTCH, RKCH5.TR5) 

Rl =R*DPHIC/(DRCHB/(1./RKCH1 +1./RKCH2) + DRCOAT/ 
K1./RKCOA1 +1./RKCOA2)) 

R2 = DRCOAT/(2.*(R-DRCOAT/2.)*DPHIB)*(1./RKCOA1 +1./RKCOA3) 
R3 = R*DPHIL/(DRCHB/(1 ./RKCH1 + 1./RKCH4) + 
1DRCOAT/(1./RKCOA1 +1./RKCOA4)) 

R4 = DRCHB/(2.*(R + DRCHB/2.)*DPHIB)*(1./RKCH1 +1./RKCH5) 
T(I,J,2) = ((T(I + 1,J,2)/R1 + T(I,J + 1,2)/R2 + T(I-1,J,2)/R3 + T{I,J-1,2) 
1/R4)/(1./R1 + 1./R2 + 1./R3 + 1 ,/R4)-T(l,J,2))*OMEGA +T(I,J,2) 

ICH(I,J) =310 
GO TO 1002 

LEFT BOUNDARY NODES 

160 IF(J.EQ.NR1)GO TO 245 
IF(J.EQ.NR2)GO TO 265 
TR1 = T(I,J,1) 

TR2 — T(I,J + 1,1) 

TR3 = T(I + 1 , J , 1 ) 

TR4 = T(I,J-1.1) 

IF(J.GT NR1.AND.J.LT.NR3)GO TO 161 

IF(J.GT.NR3)GO TO 163 

CALL METAL (MTCLO, RK1.TR1) 

CALL METAL (MTCLO, RK2.TR2) 

CALL METAL (MTCLO, RK3.TR3) 

CALL METAL (MTCLO, RK4.TR4) 

GO TO 162 

161 CALL METAL (MTCH, RK1.TR1) 

CALL METAL (MTCH, RK2.TR2) 

CALL METAL (MTCH, RK3.TR3) 

CALL METAL (MTCH, RK4.TR4) 

GO TO 162 

163 CALL METAL (MTCOAT, RK1.TR1) 

CALL METAL (MTCOAT, RK2.TR2) 

CALL METAL (MTCOAT, RK3.TR3) 

CALL METAL (MTCOAT, RK4.TR4) 

162 Rl =1./C2*(1,/RK1 +1./RK2) 
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R2 = 1./2./C1*(1./RK1 +1./RK3) 

R3 = 1./C3*(1./RK1 + 1./RK4) 

T(I,J,2) = «T(I,J + 1 ,2)/R1 + T{l + 1,J,2)/R2 + T{I,J-1,2)/R3)/ 

1(1./R1 + 1./R2 4- 1./R3)-T(l,J,2))*OMEGA-t-T(l,J,2) 

ICH(I,J) = 160 
GO TO 1002 

CHANNEL, UPPER LEFT SIDE CORNER NODE 

(ITERFACE BETWEEN COPPER AND NICKEL): ALSO CIRCUM INTERFACE 

170 R = RCO 
TR1 =T(I,J,1) 

TR2 = T(I + 1,J,1) 

TR3 = T(I,J + 1,1) 

TR4 = T{I-1,J,1) 

TR5 ~T(I,J-1,1) 

CALL METAL (MTCH, RKCH1.TR1) 

CALL METAL {MTCLO, RKCL01.TR1) 

CALL METAL (MTCLO, RKCL02.TR2) 

CALL METAL (MTCH, RKCH3.TR3) 

CALL METAL (MTCH, RKCH4.TR4) 

CALL METAL (MTCLO, RKCL04.TR4) 

CALL METAL (MTCLO, RKCL05.TR5) 

R1 =R’DPHIC/DRCL0'(1,/RKCL01 + 1./RKCL02) 

R2 = DRCHT/(R-DRCHT/2.)/DPHIL*(1./RKCH1 + 1 ./RKCH3) 

R3 = R’DPHI L/(DRCLO/(1 ./RKCLOI + 1 ./RKCL04) + DRCHT/(1 ./RKCH1 + 1 /RKCH4)) 
R4 = DRCLO/(2.*(R + DRCLO/2.)*DPHIB)*(1 ,/RKCLOI + 1 ./RKCL05) 

R5 = 2./((R*DPHIC*HC3) + (DRCHT*HC2)) 

T(I,J,2)=((T(I + 1,J,2)/R1 +T(I,J + 1,2)/R2 + T(I-1,J,2)/R3 + T(I,J-1,2)/R4 
1 +TC/R5)/(1,/R1 + 1./R2 + 1./R3 + 1./R4 + 1./R5)-T(l,J,2))*OMEGA + T(I,J,2) 

ICH(I.J) = 170 
GO TO 1002 

UPPER CHANNEL WALL 

180 R = RCO 
TR1 =T(I,J,1) 

TR3 = T(I,J-1,1) 

TR4 = T(I-1,J,1) 

CALL METAL (MTCLO, RK1.TR1) 

CALL METAL (MTCLO, RK3.TR3) 

CALL METAL (MTCLO, RK4.TR4) 

R2 = DRCLO/2./(R + DRCLO/2.)/DPHI’(1./RK1 + 1 ,/RK3) 

R3 = R*DPHI/DRCLO’(1./RK1 + 1 ,/RK4) 

R4 = 1./HC3/R/DPHI 
IF(I.EQ.NPHITOT)GO TO 185 
TR2 = T(I + 1,J,1) 

CALL METAL (MTCLO, RK2.TR2) 

R1 = R‘DPHI/DRCLO*(1./RK1 + 1./RK2) 

T(I,J,2) =((T(I + 1,J,2)/R1 +T(I,J-1,2)/R2 +T(I-1,J,2)/R3 + TC/R4)/ 

1(1. /R1 + 1./R2 + 1 ,/R3 + 1./R4)-T(l,J,2))'OMEGA + T(I,J,2) 

ICH(I,J) = 180 
GO TO 1002 

C CENTERLINE OF TOP OF CHANNEL 
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185 R2=2.*R2 
R4 = 2.*R4 

T(I,J,2) = ((T(I,J-1,2)/R2 + 1.*T(I-1,J,2)/R3 +TC/R4)/ 

1(1./R2 + 1./R3 + 1,/R4)-T(l,J,2))’OMEGA + T(I,J,2) 

ICH(I.J) = 185 
GO TO 1002 

CHANNEL SIDE WALL: ALSO A CIRCUMFERENTIAL INTERFACE 

200 TR1 — T(I,J,1) 

TR2=T(I,J-1,1) 

TR3 = T(I-1,J,1) 

TR4 = T(I,J + 1,1) 

CALL METAL (MTCH, RK1.TR1) 

CALL METAL (MTCH, RK2.TR2) 

CALL METAL (MTCH, RK3.TR3) 

CALL METAL (MTCH, RK4.TR4) 

R1 = 1 /C3L*(1,/RK1 + 1./RK2) 

R2 = 1./2./C1L*(1./RK1 + 1./RK3) 

R3 = 1,/C2L*(1./RK1 + 1./RK4) 

R4 = 1./DR/HC2 

T(I,J,2) = ((T(I,J-1.2)/R1 + T(I-1,J,2)/R2 + T(I,J + 1,2)/R3 + TC/R4)/ 

1(1 ,/R1 + 1./R2 + 1./R3 + 1./R4)-T(l,J,2))*OMEGA + T(I,J,2) 

ICH(I,J) = 200 
GO TO 1002 

WITHIN CHANNEL 
0 T(I,J,2) =TC 
ICH(I.J) =210 
GO TO 1002 

CHANNEL, LOWER LEFT SIDE CORNER NODE: 

ALSO CIRCUMFERENTIAL INTERFACE 


220 R = RCI 
TR1 =T(I,J,1) 

TR2 = T(I + 1,J,1) 

TR3 = T(I.J + 1,1) 

TR4 = T(I-1,J,1) 

TR5 = T(I,J-1,1) 

CALL METAL (MTCH, RK1.TR1) 

CALL METAL (MTCH. RK2.TR2) 

CALL METAL (MTCH, RK3.TR3) 

CALL METAL (MTCH, RK4.TR4) 

CALL METAL (MTCH, RK5.TR5) 

R1 = R*DPHIC/DRCHB*(1./RK1 + 1./RK2) 

R2 = DRCHB/(2.*(R-DRCHB/2.)*DPHIB)*(1./RK1 + 1./RK3) 

R3 = R"DPHIL/(DRCHT + DRCHB)*(1,/RK1 + 1./RK4) 

R4 = DRCHT/(R + DRCHT/2.)/DPHIL*(1,/RK1 + 1./RK5) 

R5 = 2./((R'DPHIC*HC1) + (DRCHT'HC2)) 

T(I,J,2) =((T(I + 1,J,2)/R1 +T(I,J + 1,2)/R2 + T(I-1,J,2)/R3 +T(I,J-1,2)/R4 
1 +TC/R5)/(1./R1 + 1./R2 + 1./R3 + 1./R4 + 1./R5)-T(l,J,2))*OMEGA + T(l,J,2) 
ICH(I,J) = 220 
GO TO 1002 
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LOWER CHANNEL WALL 

230 R = RCI 
TR1 =T(I,J,2) 

TR3 = T(I,J + 1,2) 

TR4 — T(I-1,J,2) 

CALL METAL (MTCH, RK1.TR1) 

CALL METAL (MTCH, RK3.TR3) 

CALL METAL (MTCH, RK4,TR4) 

R2 = DRCHB/2./(R-DRCHB/2.)/DPHI*(1./RK1 +1./RK3) 

R3 = R*DPHI/DRCHB*(1./RK1 +1./RK4) 

R4 = 1 ./HC1/R/DPHI 
IF(I EQ NPHITOT)GO TO 235 
TR2 = T(I T1.J.2) 

CALL METAL (MTCH, RK2.TR2) 

R1=R*DPHI/DRCHB*(1./RK1 + 1./RK2) 

T(I,J,2) = ((T(I + 1,J,2)/R1+T(I,J + 1,2)/R2 + T(I-1,J,2)/R3 + TC/R4)/ 
1(1. /R1 + 1./R2 + 1 ,/R3 + 1./R4)-T(l,J,2))*OMEGA +T(I,J,2) 

ICH(I,J) = 230 
GO TO 1002 

BOTTOM CENTERLINE OF CHANNEL 
5 R2 = 2 *R2 
R4 = 2.*R4 

T(I,J,2) = ((T(I,J + 1,2)/R2 + 1.*T(I-1,J,2)/R3 + TC/R4)/ 

1(1./R2 + 1./R3 + 1./R4)-T(l,J,2))*OMEGA + T(l,J,2) 

ICH(I,J) = 235 
GO TO 1002 

RIGHT BOUNDARY NODES 

250 TR1 =T(I,J,1) 

TR2 = T(I,J + 1,1) 

TR3 = T(I-1,J,1) 

TR4 = T(I,J-1,1) 

IF(J.GT.NR1.AND.J LT.NR3)GO TO 251 

IF(J GT NR3)GO TO 253 

CALL METAL (MTCLO, RK1.TR1) 

CALL METAL (MTCLO, RK2.TR2) 

CALL METAL (MTCLO, RK3.TR3) 

CALL METAL (MTCLO, RK4.TR4) 

GO TO 252 

251 CALL METAL (MTCH, RK1.TR1) 

CALL METAL (MTCH, RK2.TR2) 

CALL METAL (MTCH, RK3,TR3) 

CALL METAL (MTCH, RK4.TR4) 

GO TO 252 

253 CALL METAL (MTCOAT, RK1.TR1) 

CALL METAL (MTCOAT, RK2.TR2) 

CALL METAL (MTCOAT, RK3.TR3) 

CALL METAL (MTCOAT, RK4.TR4) 

252 R1 = 1./C2’(1./RK1 + 1./RK2) 

R2 = 1./2./C1*(1./RK1 + 1./RK3) 

R3 = 1./C3*(1./RK1 + 1./RK4) 

T(I,J,2) = ((T(I,J + 1,2)/R1 +T(I-1,J,2)/R2 + T(I,J-1,2)/R3)/ 
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1(1./R1 + 1./R2 + 1./R3)-T(l,J,2))*OMEGA + T(I,J,2) 

ICH(I,J) =250 
GO TO 1002 

INTERFACE BETWEEN COPPER AND NICKEL 

240 R = RCO 
TR1 =T(I,J,1) 

TR2 = T(I + 1,J,1) 

TR3 = T(I,J + 1,1) 

TR4 — T(I-1,J,1) 

TR5=T{I,J-1,1) 

CALL METAL (MTCH, RKCH1.TR1) 

CALL METAL (MTCLO, RKCLOI, TR1) 

CALL METAL (MTCH, RKCH2.TR2) 

CALL METAL (MTCLO, RKCL02.TR2) 

CALL METAL (MTCH, RKCH3.TR3) 

CALL METAL (MTCH, RKCH4.TR4) 

CALL METAL (MTCLO, RKCL04.TR4) 

CALL METAL (MTCLO, RKCL05.TR5) 

R1 = R* D PH l/(DRCLO/(1 /RKCLOI + 1./RKCL02) + DRCHT/(1./RKCH1 + 1./RKCH2)) 
R2 = DRCHT/2./(R-DRCHT/2 )/DPHI'(1./RKCH1 4- 1./RKCH3) 

R3 = R’DPHI/(DRCLO/(1 /RKCLOI + 1 ./RKCL04) + DRCHT/(1./RKCH1 + 1./RKCH4)) 
R4 = DRCLO/2./(R + DRCLO/2.)/DPHI*(1. /RKCLOI + 1./RKCL05) 

T(I,J,2)=((T(I + 1,J,2)/R1 + T(I,J + 1,2)/R2 + T(I-1,J,2)/R3 + T(I,J-1,2) 

1/R4)/(1./R1 + 1./R2 + 1./R3 + 1./R4)-T(l,J,2))*OMEGA + T(l,J,2) 

ICH(I,J) = 240 
GO TO 1002 

LEFT SIDE NODE (INTERFACE BETWEEN COPPER AND NIKEL) 

5 R = RCO 
TR1 =T(I,J,1) 

TR2 = T(I + 1,J,1) 

TR3 = T(I,J + 1,1) 

TR5 = T(I,J-1,1) 

CALL METAL (MTCH, RKCH1.TR1) 

CALL METAL (MTCLO, RKCLOI, TR1) 

CALL METAL (MTCH, RKCH2.TR2) 

CALL METAL (MTCLO, RKCL02.TR2) 

CALL METAL (MTCH, RKCH3.TR3) 

CALL METAL (MTCLO, RKCL05.TR5) 

R1 = R*DPHI/(DRCLO/(1 ./RKCLOI + 1 /RKCL02) + DRCHT/(1./RKCH1 + 1./RKCH2)) 
R2 = DRCHT/1./(R-DRCHT/2.)/DPHr(1./RKCH1 +1./RKCH3) 

R4 = DRCLO/1 ./(R + DRCLO/2.)/DPHI'(1. /RKCLOI + 1 ./RKCL05) 

T(I,J,2) =((1.*T(I + 1,J,2)/R1 +T(I,J + 1,2)/R2 + T(I,J-1,2)/R4) 

1/(1 ,/R1 + 1./R2 + 1./R4)-T(l,J,2))*OMEGA + T(I,J,2) 

ICH(I,J) = 245 
GO TO 1002 

INTERFACE BETWEEN TWO LAYERS WITH DIFFERENT 
RADIAL INCREMENTS (CHANNEL REGION) 

260 R = RCI 
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TR1 =T(I,J,1) 

TR2 = T{I + 

TR3 =T(I,J + 1,1) 

TR4 = T(I-1,J,1) 

TR5=T(I,J-1,1) 

CALL METAL (MTCH, RKCH1.TR1) 

CALL METAL (MTCH, RKCH2.TR2) 

CALL METAL (MTCH, RKCH3.TR3) 

CALL METAL (MTCH, RKCH4.TR4) 

CALL METAL (MTCH, RKCH5.TR5) 

R1 = R*DPHI/(DRCHT 4- DRCHB)’(1./RKCH1 + 1./RKCH2) 

R2 = DRCHB/2./(R-DRCHB/2.)/DPHI*(1./RKCH1 + 1./RKCH3) 

R3 = R*DPHI/(DRCHT + DRCHB)*(1./RKCH1 +1./RKCH4) 

R4 = DRCHT/2./(R 4- DRCHT/2.)/DPHI*(1,/RKCH1 + 1 ./RKCH5) 
T(I,J,2) = ((T(I + 1,J,2)/R1+T(I,J + 1,2)/R2 + T(I-1,J,2)/R3 + T(I,J-1,2) 
1/R4)/(1./R1 + 1 ,/R2 + 1 ,/R3 4- 1 ,/R4)-T(l,J,2))*OMEGA + T(I,J,2) 
ICH(I.J) =260 
GO TO 1002 

LEFT BOUNDARY: CHANNEL TOP/BOTTOM INTERFACE 
(DIFFERENT RADIAL INCREMENTS) 

265 R = RCI 
TR1 =T(I,J,1) 

TR2 = T(I + 1,J,1) 

TR3 = T(I,J 4-1,1) 

TR5 = T(I,J-1,1) 

CALL METAL (MTCH. RKCH1.TR1) 

CALL METAL (MTCH, RKCH2.TR2) 

CALL METAL (MTCH, RKCH3.TR3) 

CALL METAL (MTCH, RKCH5.TR5) 

R1 =R*DPHI/(DRCHT 4-DRCHB)*(1 /RKCH1 + 1./RKCH2) 

R2 = DRCHB/1 ,/(R-DRCHB/2.)/DPHI*(1 ./RKCH1 + 1 ./RKCH3) 

R4 = DRCHT/1./(R 4- DRCHT/2.)/DPHI*(1 ./RKCH1 + 1 ./RKCH5) 
T(I,J,2) = {(1.'T(I 4- 1 ,J,2)/R1 4- T(I,J 4- 1 ,2)/R2 4- T(l, J-1 ,2)/R4)/ 

1(1./R1 4-1./R24- 1./R4)-T(l,J,2))*OMEGA4T(l,J,2) 

ICH(I.J) =265 
C 

1002 RADP(I,J) = 12. *R 
IF(T(l,J,2),EQ.0.)GO TO 1003 
ERR = ABS((TT2-T(I,J,2))/T(I,J,2)) 

GO TO 1005 

1003 ERR = ABS(T(I,J,2)-TT2) 

1005 IF(ERR.GT.ERROR)IERR = IERR + 1 

1001 CONTINUE 

1000 CONTINUE 
C END OF CELL LOOP 
C WRITE(6,1006)ITER,IERR 

1006 FORMAT(10X, 'ITERATION NUMBER = ',I6,3X,'IERR = ',16) 

C DO 1007 KM = 1, NR 

C WRITE(6,2006)(T(MM,KM,2),MM = 1.NPHITOT) 

C IF(ITER1.EQ.1)WRITE(6,2003)(ICH(MM,KM),MM = 1.NPHITOT) 
C1007 CONTINUE 

2006 FORMAT(5X,11(F8.2,2X)) 



IF(IERR.GT.O)GO TO 1 
C WRITE(6,2004)ITER,ITER1 

2004 FORM AT(/10X, 'TEMPERATURE DISTRIBUTION', 5X, 

1'NUMBER OF ITERATIONS(ITER) =',I5/5X, 

2'CONDUCTIVITY ITERATIONS(ITERI) =', 15, /5X, 'RADIUS') 

C PRINT STATEMENTS TO CHECK PROCESS OF CONVERGENCE 

C DO 2000 K = 1,NR 

C WRITE(6,2001)RADP(1,K),(T(M,K,1),M = 1,NPHITOT) 

C2000 CONTINUE 
C DO 2111 K = 1,NR 

C WRITE(6,2001)RADP(1,K),(T(M,K,2),M = 1,NPHITOT) 

C2111 CONTINUE 
IEE = 0 

C TO CHECK FOR CONVERGENCE BETWEEN DIFFERENT CONDUCTIVITIES 

DO 4100 13 = 1.NPHITOT 
DO 4101 J3 = 1,NR 

EE = ABS((T(I3,J3,2)-T(I3,J3,1))/T(I3,J3,2)) 

IF(EE.GT.ERROR)IEE = IEE + 1 
4101 CONTINUE 
4100 CONTINUE 
C 

IF(ITER1.GT.50)GO TO 5 
IF(IEE.GT.O)GO TO 4 

2001 FORMAT(5X,F7.4,2X,11(F8.2,1X)) 

5 WRITE(6,3) 

IF (COOL ,EQ. 2) GOTO 4700 
WRITE(6,2005) 

WRITE (6,2007) (DPHIC*NPHIC + DPHIL*NPHIL-DPHIL*JJ,JJ = 1,NPHIL), 
1(DPHIC'NPHIC-DPHIC*II,II = 1.NPHIC) 

GOTO 4600 
4700 WRITE (6,2008) 

4600 CONTINUE 

C PRINTING THE RESULTS 

DO 2002 11 = 1, NR 

WRITE(6,2001)RADP(1,II),(T(JJ,II,1),JJ = 1.NPHITOT) 

C WRITE(6,2003)(ICH(JJ,II),JJ = 1.NPHITOT) 

2002 CONTINUE 

C DO 2222 11 = 1, NR 

C WRITE(6,2001)RADP(1,II),(T(JJ,II,2),JJ = 1.NPHITOT) 

C2222 CONTINUE 
C 

C DO 2010 MM = 1, NR 

C WRITE(6,2011)(RADP(NN,MM),NN = 1.NPHITOT) 

C2010 CONTINUE 
C 

3 FORMAT(///) 

2003 FORMAT(10X, 11(13, 2X)) 

2005 FORMAT(1H1///30X, TEMPERATURE DISTRIBUTION7/24X, 

1'RADIANS FROM CENTERLINE OF CHANNEL') 

2007 FORMAT (14X,11(F8. 5, 1X)/5X, 'RADIUS') 

2008 FORMAT (1H1///30X, 'TEMPERATURE DISTRIBUTION'//) 

2011 FORMAT(/5X,11(F8.5,2X)) 

C 
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C WRITING TEMPERATURES TO A FILE FOR CONTOURING PACKAGE 
C 

IF (IFILE EQ. 0) GOTO 5555 
WRITE (9,5100) NPHITOT,NR,NCON 
5100 FORMAT (3X,3(I3,2X)) 

WRITE (9,5200) (IFLAG(I),I = 1,7) 

5200 FORMAT (3X,7(I3,2X)) 

WRITE (9,5300) (DOT(l),l = 1,8) 

5300 FORMAT (3X,8(F7.3,2X)) 

WRITE (9,5400) IDASH(I) 

5400 FORMAT (3X,I2) 

IF (NCON ,LE. 0)NCON = -1*NCON 
DO 5500 1=1, NCON 
WRITE (9,5600) FC(I) 

5500 CONTINUE 
5600 FORMAT (3X.F7.1) 

DO 5700 II = 1.NPHITOT 
DO 5800 JJ = 1,NR 
Nl = II 
NJ=JJ 
C 

IF (T(II,JJ,1).NE.TC) GOTO 5010 
C TOP WALL OF CHANNEL 

IF (JJ.NE.NR1 + 1) GOTO 5020 
NJ = NR1 
Nl = NPHIL 
GOTO 5010 
5020 CONTINUE 

C BOTTOM WALL OF CHANNEL 

IF (JJ.NE.NR2-1) GOTO 5030 
NJ = NR2 
Nl = NPHIL 
GOTO 5010 
5030 CONTINUE 

C SIDE WALL AND INTERIOR OF CHANNEL 
Nl = NPHIL 
C 

5010 CONTINUE 

C 

IF (NI.GT. NPHIL) GOTO 5040 
PSI =(NI-1)*DPHIL 
5040 CONTINUE 

IF (Nl. LT. NPHIL + 1) GOTO 5050 
PSI = (N PHIL-1 )*DPHIL + DPHIC*(NI-N PHIL) 

5050 CONTINUE 

C 

Zl = RADP(NI,NJ)'SIN(PSI) 

ZJ = RADP(NI,NJ)'COS(PSI) 

C 

WRITE (9,5900) ZI,ZJ,T(NI,NJ,1) 

C 

C 

C 

5900 FORMAT (3X,2(F7.4,2X),F8.2) 
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5800 CONTINUE 
5700 CONTINUE 
C 

5555 CONTINUE 

CALCULATE HEAT TRANSFER RATES 
THESE ARE BULK RATES FOR EACH CHANNEL USING FEET 
QO = 0. 

QC = 0. 

QC1 =0. 

QC2 =0. 

QC3 = 0. 

Ql =0. 

HOAVE =0. 

DO 3000 1 = 1, NPHITOT 
HOAVE = HOAVE + HO(l) 

DO 3001 J = 1,NR 
IF (I.LT.NPHIL)DPHI = DPHIL 
IF (I.GT.NPHIL)DPHI = DPHIC 
IF(ICH(I,J).EQ.320)GO TO 3700 
IF(ICH{I,J).EQ.300)GO TO 3750 
IF(ICH(I,J).EQ.230)GO TO 3100 
IF(ICH(I,J).EQ.180)GO TO 3150 
IF(ICH(l,J).EQ.235)GO TO 3200 
IF(ICH(l,J).EQ.185)GO TO 3250 
IF(ICH(I,J) EQ.200)GO TO 3300 
IF(ICH(I,J).EQ.220)GO TO 3400 
IF(ICH(I,J),EQ.170)GO TO 3450 
IF(ICH(I,J),EQ,120)GO TO 3500 
IF(ICH(l,J),EQ.100.OR.ICH(l,J).EQ.110)GO TO 3550 
IF(ICH(I,J) EQ.150)GO TO 3600 
IF(ICH(l,J).EQ.130.OR.ICH(l,J).EQ.140)GO TO 3650 
GO TO 3001 

3100 QC1 =QC1 +2.*RCI*DPHI*HC1*(T(I,J,2)-TC)*1. 

GO TO 3001 

3150 QC3 = QC3 + 2.*RCO*DPHI*HC3*(T(l,J,2)-TC)*1. 

GO TO 3001 

3200 QC1 =QC1 + RCI*DPHIC*HC1*(T(I,J,2)-TC)*1. 

GO TO 3001 

3250 QC3 = QC3 + RCO'DPHIC'HC3'(T(l,J,2)-TC)*1. 

GO TO 3001 

3300 QC2 = QC2 + 2 •DRCHT*HC2*(T(I,J,2)-TC)*1. 

GO TO 3001 

3400 QC1 =QC1 +(RCPDPHIC*HC1)*(T(I,J,2)-TC)*1. 

QC2 = QC2 + (DRCHT”HC2)*(T(I,J,2)-TC)'1. 

GO TO 3001 

3450 QC3 = QC3 + (RCO*DPHIC*HC3)*(T(l,J,2)-TC)*1 . 

QC2 = QC2 + (DRCHT’HC2)*(T(I,J,2)-TC)*1. 

GO TO 3001 

3500 QO = QO + 2.*RO’DPHrHO(l)*(T{l,J,2)-TO)*1. 

GO TO 3001 

3550 QO = QO + RO*DPHI*HO(l)*(T(I.J,2)-TO)*1 . 

GO TO 3001 

3600 Ql = Ql + 2.*(RI-TCOAT)*DPHI*HI*(T(l,J,2)-TI)*1 . 
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GO TO 3001 

3650 QI = QI +(RI-TCOAT)*DPHrHr(T(l,J,2)-TI)*1. 

GO TO 3001 

3700 QI = QI + 2.*(RI-TCOAT)'DPHIB*HI*(T(l,J,2)-TI)*1. 

GO TO 3001 

3750 QO — QO + 2.*RO*DPHIB*HO(l)*(T(l,J,2)-TO)'1. 

3001 CONTINUE 
3000 CONTINUE 

QC = QC1 + QC2 + QC3 

Q = QC + QO + QI 

HOAVE = HOAVE/NPHITOT 

WRITE(6,3002)QC,QC1,QC2,QC3,QO,QI,Q 

3002 FORMAT(1H1/10X,'HEAT RATE PER CHANNEL IN BTU/S PER FT LENGTH' 
1/10X, 'COOLING CHANNEL = ',F13.2/10X, 'INNER WALL = ',F13.2,2X, 

2'SIDE WALL = ',F13.2,2X,'OUTER WALL = '.F13.2/10X, 

3'HEAT TRANSFER TO OUTSIDE = '.F13.2 
4/10X/HEAT TRANSFER FROM HOT GAS WALL = ',F13.2 
5/10X/SUM OF ALL BULK HEAT TRANSFERS = ',F13.2) 

WRITE(6,3) 

C CALCULATING THE AVERAGE FLUXES FOR DIFFERENT WALL SECTIONS 
C THESE ARE IN BTU/IN -2-S) 

QO = QO*NCC/(2.*3.1415‘RO)/144. 

QI = QI'NCC/(2. , 3.1415*(RI-TCOAT))/144. 

IF (COOL EQ. 2) GOTO 3010 
QC1 = QC1/(2."DPHIC*NPHIC’RCI)/144. 

QC2 = QC2/(2.*(RCO-RCI))/144. 

QC3 = QC3/(2.’DPHIC*NPHIC*RCO)/144. 

QC = QC/(2.*(DPHIC*NPHIC*(RCI + RCO) + (RCO-RCI)))/144. 

3010 CONTINUE 
C 

WRITE(6,3005)QI,QC,QC1,QC2,QC3,QO 
3005 FORMAT(////10X,'HEAT FLUXES AT THIS STATION IN BTU/IN _, 2-S' 
1/10X/GAS WALL = ',F13.2/10X,' AVER AGE COOLANT CHANNEL = ',F13.2 
2/10X, 'INNER WALL = ',F13.2,2X,'SIDE WALL = ',F13.2,2X, 

3'OUTER WALL = '.F13.2/10X, 'CLOSE-OUT = '.F13.2) 

C 

WRITE(6,3) 

WRITE(6,4500)HOAVE,ITER1,ITER2 

4500 FORM AT(///10X, 'OUTSIDE HEAT TRANSFER COEFFICIENT (BTU/FT -2-) = ', 
1F10.5/10X, 'CONDUCTIVITY ITERATIONS(ITERI) =',I5/10X, 

2'TOTAL NUMBER OF ITERATIONS = ',15) 

RETURN 

END 

C 

SUBROUTINE COPPER(RKC.TC) 

SUBROUTINE FOR EVALUATION OF COPPER THERMAL CONDUCTIVITY 

NINT NUMBER OF INTERVALS 

RK CONDUCTIVITY 

T TEMPERATURE 

DIMENSION RK(10),T(10) 

NINT = 2 



ooooooo o ooonnoo 


RK(1) = 53.*10.**(-4)*12. 

T(1) = 400. 

RK(2) = 52.*10.**(-4)*12. 

T(2) = 500. 

RK(3) =46.6*10.**(-4)*12. 

T{3) = 1650. 

CALL INTER(RK,T,NINT,TC,RKC) 

RETURN 

END 

C 

SUBROUTINE NICKEL(RKC.TC) 

SUBROUTINE FOR EVALUATION OF NICKEL THERMAL CONDUCTIVITY 

NINT NUMBER OF INTERVALS 

RK CONDUCTIVITY 

T TEMPERATURE 

DIMENSION RK(10),T(10) 

NINT = 5 

RK(1)=26.*10.**(-4)*12. 

T( 1 ) = 1 30. 

RK(2) = 15.*10.**(-4)*12. 

T(2) = 300. 

RK(3) =9.5*10.**(-4)*12. 

T(3) =500. 

RK(4) = 7.2*10.**(-4)*12. 

T(4) = 800. 

RK(5) =6.2*10.**(-4)*12. 

T{5) = 1200. 

RK(6) = 7.2*10.*’(-4)*12. 

T(6) = 1800. 

CALL INTER(RK,T,NINT,TC,RKC) 

RETURN 
END 

SUBROUTINE SOOT(RKC.TC) 

SUBROUTINE FOR EVALUATION OF SOOT CONDUCTIVITY 

NINT NUMBER OF INTERVALS 

RK CONDUCTIVITY 

T TEMPERATURE 

DIMENSION RK(20),T(20) 

NINT = 1 

RK(1) = 7.’10.**(-6)*12. 

T(1) = 100. 

RK(2) = 7.*10.’*(-6)*12. 

T(2) = 2000. 

CALL INTER(RK,T,NINT,TC,RKC) 

RETURN 
END 
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ooooooo n nnnnnnn oo ooooooo 


SUBROUTINE NARLOYZ(RKC.TC) 

SUBROUTINE FOR EVALUATION OF NARLOY-Z CONDUCTIVITY 

NINT NUMBER OF INTERVALS 

RK CONDUCTIVITY 

T TEMPERATURE 

DIMENSION RK(20),T(20) 

NINT = 1 

RK(1) = 4.0375*10.**{-3.)*12. 

T(1) = 500. 

RK(2) = 4.305*10.**(-3.)*12. 

T(2) = 1000. 

CALL INTER(RK,T,NINT,TC,RKC) 

RETURN 

END 


SUBROUTINE RSR995(RKC,TC) 

SUBROUTINE FOR EVALUATION OF RSR 995-AE CONDUCTIVITY 

NINT NUMBER OF INTERVALS 

RK CONDUCTIVITY 

T TEMPERATURE 

DIMENSION RK(20),T(20) 

NINT = 4 

RK(1) =4.95’10.**(-3.)*12. 

T(1) = 100. 

RK(2) =4.95*10.‘*(-3.)*12. 

T{2) = 300. 

RK(3) =4.86*10.**(-3.)*12. 

T(3) = 500. 

RK(4) =4.699 *10. **(-3.) '12. 

T(4) = 1000. 

RK(5) = 4.37*10.**(-3.)*12 
T(5) = 1500. 

CALL INTER(RK,T,NINT,TC,RKC) 

RETURN 

END 

SUBROUTINE COLUMB(RKC.TC) 

SUBROUTINE FOR EVALUATION OF COLUMBIUM CONDUCTIVITY 

NINT NUMBER OF INTERVALS 

RK CONDUCTIVITY 

T TEMPERATURE 

DIMENSION RK(20),T(20) 

NINT = 1 

RK(1)=4.5*10.**(-4)*12. 
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non no 


T(1) = 100. 

RK(2) = 4.5*10.**(-4)*12. 

T(2) = 2000. 

CALL INTER(RK,T,NINT,TC,RKC) 

RETURN 

END 


SUBROUTINE INTER(RK,T,NINT,TC,RKC) 

SUBROUTINE FOR INTERPOLATION 

DIMENSION RK(10),T(10) 

DO 10 I = 1.NINT 
M = I 

IF(TC.GE.T(I).AND.TC.LE.T(I -F 1))GO TO 20 
10 CONTINUE 
GO TO 30 

20 RKC = RK(M) + (RK(M + 1)-RK(M))*(TC-T(M))/(T(M + 1)-T(M)) 

GO TO 40 

30 IF(TC.LT.T(1))RKC = RK(1) + (RK(2)-RK(1))*(TC-T(1)) 

1 / (T (2)-T ( 1 )) 

IF(TC.GT T(NINT + 1))RKC = RK(NINT) + (RK(NINT + l)-RK(NINT))* 
1(TC-T(NINT))/(T(NINT + l)-T(NINT)) 

40 RETURN 
END 
C 

SUBROUTINE METAL (M,RK,T) 

C THIS SUBROUTINE ALLOWS THE TRIAL OF MANY DIFFERENT MATERIALS FOR 
C THE DIFFERENT REGIONS. TO ADD A NEW MATERIAL, ADD A NEW NUMBER IN 

C THE GOTO LIST, THE SUBROUTINE TO FIGURE CONDUCTIVITY AND THE CALL 

C TO THAT ROUTINE. MAKE SURE TO CHANGE THE MATERIAL NUMBER. 

C 

GOTO (1,2, 3, 4, 5, 6), M 

1 CALL COPPER(RK, T) 

RETURN 

2 CALL NICKEL(RK.T) 

RETURN 

3 CALL SOOT(RK.T) 

RETURN 

4 CALL NARLOYZ(RK.T) 

RETURN 

5 CALL RSR995(RK,T) 

RETURN 

6 CALL COLUMB(RK.T) 

RETURN 

END 

/EOF 



Sample Output 


CASE SRB 


NPHIL = 4 NPHIC = 3 

CW = 0. 031415 DIFPHIL = 0.025733 DIFPHIC = 0.005682 
DPHIL = 0. 008578 DPHIC = 0.001894 DPHIB = 0.005236 


NRCLO = 5 NRCHT = 5 NRCHB = 4 NRCOAT = 0 

DRCLO =0.000417 
DRCHT = 0.002750 
DRCHB =0.000417 
DRCOAT = 0.000000 
DPHIC = 0.001894 
DPHIL = 0.008578 
NR = 14 

NR1 = 5 NR2 = 10 NR3 = 14 


TEMPERATURE DISTRIBUTION 
RADIANS FROM CENTERLINE OF CHANNEL 



0.03141 0.02284 

0.01426 

0.00568 

0.00379 

0.00189 

0.00000 

1.9450 

113.48 

113.42 

113.26 

113.03 

113.00 

112.98 

112.97 

1.9400 

113.48 

113 43 

113.26 

113.02 

112.99 

112.96 

112.95 

1.9350 

113.50 

113.44 

113.27 

113.00 

112.95 

112.92 

112.91 

1.9300 

113.52 

113.46 

113.28 

112.98 

112 89 

112.83 

112.81 

1 9250 

113.56 

113.49 

113.31 

112.98 

112.77 

112.67 

112.64 

1.8920 

113.94 

113.88 

113.70 

113.39 

90.00 

90.00 

90.00 

1.8590 

114.83 

114.77 

114.59 

114.28 

90.00 

90.00 

90.00 

1.8260 

116.26 

116.20 

116.02 

115.71 

90.00 

90.00 

90.00 

1.7930 

118.27 

118.22 

118.05 

117.75 

90.00 

90.00 

90.00 

1.7600 

120.81 

120.80 

120.76 

120.71 

121.00 

121,13 

121.17 

1.7550 

121.21 

121.20 

121.18 

121.21 

121.31 

121.37 

121.39 

1.7500 

121.61 

121.61 

121.60 

121.65 

121.68 

121.71 

121.72 

1.7450 

122.02 

122.02 

122.02 

122.06 

122.08 

122.09 

122.10 

1.7400 

122.43 

122 42 

122.43 

122.47 

122.48 

122.49 

122.50 


HEAT RATE PER CHANNEL IN BTU/S PER FT LENGTH 
COOLING CHANNEL= 0.60 

INNER WALL = 0.04 SIDE WALL = 0.53 OUTER WALL = 

HEAT TRANSFER TO OUTSIDE = 0.00 

HEAT TRANSFER FROM HOT GAS WALL= -0.60 


0.03 
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SUM OF ALL BULK HEAT TRANSFERS = 


0.00 


HEAT FLUXES AT THIS STATION IN BTU/IN-2-S 

GAS WALL = -0.46 

AVERAGE COOLANT CHANNEL = 0.13 

INNER WALL= 0.15 SIDE WALL= 0.13 OUTER WALL = 
CLOSE-OUT = 0.00 


OUTSIDE HEAT TRANSFER COEFFICIENT (BTU/FT - 2-) = 0.00048 

CONDUCTIVITY ITERATIONS(ITERI) = 4 


0.12 
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FIGURE 1. - A ROCKET (SPACECRAFT) THRUST CHAMBER. 


49 



NICKEL 

I 



RADI/ 
INCRE 
IN Cl 
(DRCt 


RADIAL 
TOP CHA 


NG 

IEL 


RADIAL 
TOM Cl 



0 400 800 1200 1600 2000 2400 2800 5200 3600 

TEMPERATURE, °R 

FIGURE 4. - THERMAL CONDUCTIVITY OF NICKEL (FROM REF. 3). 
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FIGURE 5. - THERMAL CONDUCTIVITY OF COPPER (FROM REF. 3). 



0 400 800 1200 1600 2000 2400 2800 3200 3600 
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FIGURE 6. - THERMAL CONDUCTIVITY OF COLUMBIUM (FROM REF. 3). 
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